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Abstract. Causality and stability in relativistic dissipative hydrodynamics 
are important conceptual issues. We argue that causality is not restricted to 
hyperbolic set of differential equations. E.g. heat conduction equation can 
be causal considering the physical validity of the theory. Furthermore we 
propose a new concept of relativistic internal energy that clearly separates 
the dissipative and non-dissipative effects. We prove that with this choice we 
remove all known instabilities of the linear response approximation of viscous 
and heat conducting relativistic fluids. In this paper the Eckart choice of the 
velocity field is applied. 



1. Introduction 

The theories of dissipative fluids are different in the relativistic and nonrelativis- 
tic spacetime. The simplest first order, parabolic theory of nonrelativistic fluids, the 
system of Navier-Stokes and Fourier equations, is tested and justified by countless 
applications of the everyday engineering practice. The second order generalization 
of the theory introduces the fluxes of the extensives as independent variables. In 
this way the validity and applicability of the hydrodynamic and heat conduction 
equations is extended providing a hyperbolic system [HE]- Although the well known 
special relativistic generalization of the first order Navier-Stokes-Fourier system is 
straightforward [3J, it has some unacceptable features. It is acausal and unstable. 
Therefore it is generally accepted that only the extended, second order theories are 
viable. Nowadays, heavy ion collision experiments give a unique opportunity to 
check the different suggestions and the interest in dissipative relativistic theories is 
renewed [H [5] , [61 H E El HOl HH E] . 

The first order theories are based on the local equilibrium hypothesis, where the 
independent variables are the same as in equilibrium, but in second (and higher) 
order theories the fluxes of the local equilibrium theory appear as independent 
variables. An other usual property is that in case second order relativistic theories 
the entropy vector is quadratic in the fluxes, containing terms like q a q a u tJ ', IF 1 "^ 
etc., characterizing the deviation from the local equilibrium. As one can introduce 
general dynamic variables beyond the fluxes (see e.g. [13]), the above classification 
is not general. 

Differential equations of first order theories are parabolic, therefore they are con- 
sidered generally as acausal. The differential equations of the second order theories, 
that are constructed according to the Second Law, are mostly hyperbolic, therefore 
they are considered generally as causal. However, more careful considerations show, 
that the relation between parabolicity and causality is not so straightforward and 
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requires some attention both from a mathematical and from a physical point of 
view. 

The homogeneous equilibrium in first order theories is generally considered as 
unstable. The homogeneous equilibrium of second order theories is generally con- 
sidered as stable. The stability considerations are referring to the linear stability 
calculations of Hiscock and Lindblom [21 [T5] . 

However, like these second order theories are extending the validity of the first 
order ones, their physical content is included in the corresponding more general sec- 
ond order theory. The more involved second order theories do not cure necessarily 
the instabilities of the first order theories. It is shown by Geroch and Lindblom 
that physical fluid states in these theories relax to the solutions of the underlying 
first order theory [16j E] • 

The missing of a simple and stable relativistic generalization of the Navier- 
Stokes-Fourier theory resulted in several attempts to improve the properties of first 
order theories. Garcia-Colm and Sandoval- Villalbazo suggested a separation of an 
internal energy balance from the balance of the energy-momentum, similarly to 
nonrelativistic theories [18 . However, with an additional independent energy bal- 
ance the energy-momentum tensor would not embrace the whole energy content of 
the matter [TI5] . Other authors suggest a suitable definition of the four- velocity field 
[2"0] . None of the previous suggestions investigate the stability of the corresponding 
equations. 

In the following section we argue that the speed of the propagation of signals can 
be finite in parabolic theories, too, if their physical validity is considered. Therefore 
first order relativistic hydrodynamic theories cannot be excluded by referring to 
causality. Moreover, in the light of the above mentioned observation of Geroch 
and Lindblom it is even more important to find a viable relativistic generalization 
of the Navier-Stokes-Fourier equations. A necessary condition for a causality in a 
weaker sense [21] is the stability of the homogeneous solutions of the corresponding 
differential equations. Based on this observation in the followings we outline a 
new approach to relativistic fluids. We suggest a separation of the dissipative and 
non-dissipative parts of the energy momentum distinguishing between the total 
energy density and the internal energy density of the matter. The later is the 
absolute value of the projected energy flux four-vector, this way it incorporates the 
momentum density as well. Since in the corresponding thermodynamic frame the 
entropy density depends also on the energy flux besides the energy density, but does 
not on the pressure, our suggestion can be classified between the first order, local 
equilibrium one, and of the extended, second order theories. In the final section we 
demonstrate the linear stability of the homogeneous equilibrium of viscous, heat 
conducting fluids with the Eckart choice of the velocity field. 

2. Remarks on causality 

The common argument against the use of parabolic differential equations in 
physics is that some of their typical solutions show signal propagation with infi- 
nite speed. More sophisticated arguments require a well posedness of the related 
mathematical problems that can be guaranteed by hyperbolicity. The character- 
istic surface of the simplest relativistic heat conduction equation with constant 
coefficient is a spacelike hypersurface according to a comoving observer. Moreover, 
the characteristic surfaces are invariant to the transformation of the equation (in 
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particular to a Lorentz boost of the reference frame) : The parabolicity or the hyper- 
bolicity of the equation does not change by changing the observer. Therefore, as it 
was argued by Kostadt and Liu [22] , by simple mathematics initial value problems 
of parabolic differential equations can be well posed, provided that initial data are 
given on the characteristic surface of the equations. From a physical point of view 
this is a natural requirement. 

On the other hand, the characteristic hypersurfaces are those that determine the 
speed of propagation of simple solutions (the domain of influence) of a hyperbolic 
differential equation, too. Therefore, speed of the signal propagation for a hyper- 
bolic differential equation is in general not infinite, but can be higher than the 
speed of light. The actual speed depends on the parameter values in the equation. 
This statement appears as a trivial fact in case of wave propagation equation. In 
one space dimension, considering a comoving observer, with respect to a constant 
velocity field u a we get the following form 

(1) d tt 9 - c 2 w d xx 9 = 0. 

Here 9 is the corresponding scalar physical quantity, c w is the wave propaga- 
tion speed. The solution of the characteristic differential equation of (JTJ) gives the 
equation 9(x, t) = x ± c w t = const, for the two characteristic lines of the equation. 
Applying a Lorentz transformation x = j(x — vt) and t = j(t — vx/c 2 ) we can get 
the transformed form of the characteristic lines as: 

6(x, t) = (1 ± ^)x + (v± c w )t. 

Therefore the transformed characteristic speed is 

_ _ v±c w 
Cw ~~ -y + CzuV ' 

c 2 

We can get the same result with the Lorentz transformation of the equation, 
too. The above expression shows that the propagation speed of waves can be faster 
than the speed of the light only if c w > c, as we have expected. 

In case of a set of nonlinear differential equations the calculation of the charac- 
teristic wave speeds can be more involved and even the proof of the hyperbolicity 
of the corresponding set of equations is not trivial. In general the value of the 
speed will depend on parameters in the set of equations and the relativistic, covari- 
ant form combined with hyperbolicity do not warrant a propagation speed smaller 
than the speed of light. 

On the other hand the theoretically infinite speeds in parabolic equations are 
usually not observable, because their effect is out of the physical validity range of 
the theory. The atomistic structure of the matter restricts the validity of continuum 
descriptions as it was pointed out by Weynmann [231122]. With the help of the mean 
free path and the collision time one can give simple estimates on the propagation 
speed of measurable signals. 

We demonstrate this property on the example of the Fourier heat conduction 
equation 



(2) 



8 t 9 - \d xx 9 = 0. 
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Figure 1. Mean free path limited signal propagation according 
to the Fourier heat conduction equation. £ = .2, A = 1 and t — 
0.2,0.4,0.6,0.8 



The hydrodynamic range of validity requires that 6 must not vary too rapidly 
over a mean free path £ 

1 d9_ 



(3) 



1 

K< 1 



Assuming a sharp initial condition the solution of the heat conduction equation 
can be written as 

0(x, t) — p 

This is a typical acausal solution of the Fourier equation. However, substituting 
the above solution into the condition ^ we get a limit of the propagation speed of 
the continuum signals as 

fA\ x ^ 

(4) - << Vu m oc -. 

Therefore, instead of the infinite tail of the solution in the reality we have an 
extending range, cf. (P). 

In case of heat conduction in water at room temperature we can easily give 
an estimate as v max ~ ^7 — 14m/s, where A is the Fourier heat conduction 
coefficient, p is the density and c v is the specific heat of water. As heat conduction 
is disputed phenomena in quark-gluon plasma we cannot give a reliable estimation 
here. 

Independently of the previous estimation Fichera suggested that the speed of 
the signal propagation is restricted by observability of the given physical quantity. 
Observability can be related to the sensitivity of the measurement but also to 
fluctuations and the particular structure of the matter. That can give an other 
bound to the speed of the signal propagation J24J [25l [21]. In our case we may 
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FIGURE 2. Signal propagation according to the observability 
threshold of the Fourier heat conduction equation. 9 max = 0.1, 
A = 1 and t = 0.2,0.4,0.6,0.8 



assume that we cannot observe \0\ below a given value \6\ < 9 max . Then the 
propagation speed becomes finite, nevertheless it is not constant as one can inspect 
in figure |[2}. 

Summarizing our arguments have seen that hyperbolicity of the equations can 
lead to well posed problems and gives finite propagation speed but does not warrant 
that the propagation speed is less than the speed of light. On the other hand 
one can formulate well posed Cauchy problems related to parabolic equations, too. 
Moreover, the physical validity of a continuum theory can warrant slow propagation 
speeds in several different ways. Therefore, we may conclude that parabolic and 
mixed systems of continuum differential equations (as Fourier heat conduction or 
any first order continuum hydrodynamics) could be useful models in relativistic 
theories. Those theories cannot be excluded by causality arguments. 

All the previous estimates are connected to some definite properties of the so- 
lutions of the simple heat conduction equation. If the exponential damping of 
the solutions cannot be guaranteed, then causality issues can become important. 
Therefore the stability of the homogeneous equilibrium is not only an evident phys- 
ical requirement but also a necessary condition for the causality of any first order 
dissipative relativistic hydrodynamics. 

3. Balances of particle number, energy and momentum of 

relativistic fluids 

For the metric (Lorentz form) we use the g^ v — diag(— 1, 1, 1, 1) convention 
and we use a unit speed of light c = 1, therefore for a four- velocity u a we have 
u a u a — —1. A ^ = g a p + u a Ufj denotes the w-orthogonal projection. This metric 
convention will be convenient in the stability calculations. 
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In the following we fix the velocity field to the particle number flow, according 
to Eckart. Therefore the particle number flow is timelikc by definition and can be 
expressed by the local rest frame quantities as 

(5) N a = nu a . 

Here n = —u a N a is the particle density in comoving frame. 
The particle number conservation is described by 

(6) d a N a = n + nd a u a = 0, 

where h = 4p = u a d a n denotes the derivative of n with respect to the proper time 

T. 

The energy-momentum density tensor is given with the help of the rest-frame 
quantities as 

(7) T a/3 = eu a vP + u a q^ + vPq a + P a/3 , 

where e = u a upT a P is the density of the energy, q^ — — u Q A^ 7 T a7 is the energy 
flux or heat flux, q a = -u^A^T^ is the momentum density and P af} = A"A^T W 
is the pressure (stress) tensor. The momentum density, the energy flux and the 
pressure are spacelike in the comoving frame, therefore u a q a — and u a P a @ — 
U pP a P — Ql 3 . Let us emphasize that the ((TJ) form of the energy- momentum tensor 
is completely general, it is just expressed by the local rest frame quantities. The 
energy-momentum tensor is symmetric, because we assume that the internal spin 
of the material is zero. In this case the heat flux and the momentum density are 
equal. However, the difference in their physical meaning is a key element of our 
train of thoughts. Heat is related to dissipation of energy but momentum density 
is not, therefore this difference should appear in the corresponding thermodynamic 
framework. 

Now the conservation of energy-momentum c^T"' 3 = is expanded to 

(8) d p T al3 = eu a + eu^dfiu? + eu a + u a d q f) + q^dpu* + <f + q a d f3 u 13 + d P af} . 
Its timelike part in the local rest frame gives the balance of the energy e 

(9) - u a d T a P = e + ed a u a + d a q a + q a ii a + P a0 d p u a = 0. 

The spacelikc part in the local rest frame describes the balance of the momentum 

(10) A Q 7 ^T^ = eu a + q a dfiv? + q^dpu + A%q^ + A«dpP^ = 0. 

4. Thermodynamics 

The entropy density and flux can also be combined into a four-vector, using local 
rest frame quantities: 

(11) S a = su a + J a , 

where s = —u a S a is the entropy density and J a — S a — u a s = A a pS l3 is the 
entropy flux. The entropy flux is u-spacelike, therefore u a J a — 0. Now the Second 
Law of thermodynamics is translated to the following inequality 

(12) d a S a = s + sd a u a + d a J a > 

Relativistic thermodynamic theories assume that the entropy is a function of 
the local rest frame quantities, because the thermodynamic relations reflect general 
properties of local material interactions. The most important assumption is that 
the entropy is a function of the local rest frame energy density, the time-timelike 
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component of the energy momentum tensor according to the velocity field of the 
material 3, 26J. Definitely the thermodynamics cannot be related to an external 
observer, therefore the dependence on the relative kinetic energy is excluded. This 
interpretation of e in ||7J) is supported by the form of the energy balance ([9]) , where 
the last term is analogous to the corresponding internal energy source (dissipated 
power) of the nonrelativistic theories. 

In nonrelativistic fluids the internal energy is the difference of the conserved 
total energy and the kinetic energy of the material. However, also in nonrelativistic 
theories the constitutive relations must be objective in the sense that they cannot 
depend on an external observer, the thermodynamic framework should produce 
frame independent material equations. (This apparent contradiction of classical 
physics is eliminated by different sophisticated methods and lead to such important 
concepts as the configurational forces or/and virtual power [27| . l28 l 129]). However, 
without distinguishing the energy related to the flow of the material from the total 
energy one mixes the dissipative and nondissipative effects. The wrong separation 
leads to generic instabilities of the corresponding theory. 

Our candidate of the relativistic internal energy is related to the energy vector 
defined by E a = —u^T a ^ = eu a + q a . The energy vector embraces both the total 
rest frame energy density and the rest frame momentum. Therefore its absolute 
value ||i?|| = \/—E a E a = \J\e? — q a q a \ seems to be a reasonable choice of the 
scalar internal energjQ. Its series expansion, when the energy density is larger than 
the momentum density, shows strong analogies to the corresponding nonrelativistic 
definition 



2 



\\E\\=VW=q^\*e-^ + ... 

Thermodynamic calculations based on the Liu procedure support this assump- 
tion [30] . Let us emphasize, that our candidate of internal energy is not related 
to any external reference frame, only to the velocity field of the material. In a 
Landau-Lifschitz frame the energy vector is timelike. 

Assuming that the entropy density is the function of the internal energy and the 
particle number density s(e,g Q ,n) = s(-^/|e 2 — q a q a \,n) leads to a modified form 
of the thermodynamic Gibbs relation and the potential relation for the densities as 
follows 

q a q 2 
(13) de dq a = Tds + /idn, and e = Ts — p + fin. 

e e 

Here T is the temperature, p is the pressure and fi is the chemical potential. 
Equivalently the Gibbs relation gives the derivatives of the entropy density as fol- 
lows 



ds 
de 



(g a ,n) 



I i-1 

T" dq" 



(e,n) 



<la_ 9S_ 

eT' dn 



T' 



For the entropy flux we assume the classical form 

(14) J a = jr- 

Now we substitute the energy balance (|9|) and the particle number balance ([6]) 
into the entropy balance (|12p and we arrive at the following entropy production 



Note that q a is spacelike, therefore the quantity under the sign is non-negative. 
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formula: 

d a S a = s(e,q a ,n) + sd a u a + d a J a 

ds ds ds a q a 

= -^- e +^-zQ + -^-n + sd a u + d a — 
de dq a on T 

= -^d a u a + d a q a + q a u a + P af3 d p u a ) - ^q a + sd a U a 
1 1 e 

+ t.nd a u a + q a d a ± + ^d a q a 

(15) = -^(P a e-(-e + S T + vn)A^)d a u p + q«(d a ±-^-^)>0. 

According to this quadratic expression and the potential relation in (|13|) the 

viscous pressure is given by 

We may introduce the conventional decomposition of the pressure 

n 11 " = (p + n)A Q/3 + (n Q/3 ). 

Here II = ±11% -p = - £ and (IT*' 3 ) = A£A£(|(IF" / + W) - \A^ V A 7(5 IF 5 
is symmetric an traceless. 

Therefore the (|15[) entropy production can be written as 

(16) - ^)d aUp - 1 (n- + f (d a ± - ^ - > o 

In isotropic continua the above entropy production results in the following con- 
stitutive functions assuming a linear relationship between thermodynamic fluxes 
and forces 

(17) q a = -\l 5 A^(a r T + Ti h + ^ 

(18) (U afi ) = -2n{d a u p ), 

2 

(19) II-^L = -n v 8 a u a . 

e 

Let us recognize the additional term to the bulk viscous pressure. For the whole 
viscous stress we get 

(20) n Q/3 = -2rj{d a u ) - 77„A Q/3 a 7 u 7 . 

(|17p and (|2"0|) are the relativistic generalizations of the Fourier law of heat con- 
duction and the Newtonian viscous pressure function. The shear and bulk viscosity 
coefficients, n and r/ v and the heat conduction coefficient A = AT~ 2 are non nega- 
tive, according to the inequality of the entropy production l|17p. We may introduce 
a relaxation time r = A/e in (117[) . as usual in second order theories. 

The equations ^ and §§§ are the evolution equations of a relativistic heat con- 
ducting ideal fluid, together with the constitutive function (|2"D|) and the relaxation 
type equation (|T7f . As special cases we can get the relativistic Navier-Stokes equa- 
tion substituting ([20| into (flO|) and assuming q a = 0, or the equations of rela- 
tivistic heat conduction solving together (fl"7|) and © assuming that IT" 3 = and 
u a = const.. 
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5. Linear stability 

In this section we investigate the linear stability of the homogeneous equilibrium 
of the equations ©, © and (fTU]) together with the constitutive relations ([20)1 and 
([T"7|) . Similar calculations are given by Hiscock and Lindblom both for Eckart fluids 
[14] and Israel-Stewart fluids [15]. 



5.1. Equilibrium. The equilibrium of the above set of equations is defined by 
vanishing proper time derivatives and by zero entropy production with vanishing 
thermodynamic fluxes 

(21) n Q/3 = and q a = 0. 

Therefore according to the balances and the constitutive functions the equilib- 
rium of the fluid is determined by 

(22) n = const, e = const. => T = const., /j, — const., p — const., 

(23) d a U a = 0, d a Ufj + d(3U a = 0. 

In addition to the above conditions we require a homogeneous equilibrium ve- 
locity field 

(24) u a = const. 



5.2. Linearization. We denote the equilibrium fields by zero lower index and the 
perturbed fields by 5 as Q — Q + 5Q. Here Q stands for n, e, u a , q a , and IT a/3 . 
The linearized equations d6j) , (J9j) , (flO]) , (flT]) , ([20]) around the equilibrium given by 
(ED- El)- (El become 



(25) = Sn + nd a Su a , 

(26) = 5e+ (e+p)daSu a + d a 6q a , 

(27) = (e + p)5u a + A aP d p Sp + Sq a + A^dpSW 13 , 

(28) = Sq a + AA" 7 (d~ t 5T + T6u 7 + ^Sq 7 

(29) = 5Il afl +ij v d 1 Su~<A al3 +r 1 A a ~<A f3 > 1 {d J 5u ll + d ll Su 1 ). 

Here fj v — (rj v — | rf). The perturbation variables satisfy the following properties 
inherited from the linearization of the original ones 

= U a Sq a = U a 5u a = U a 5U a l3 = SH a /3 — SH/3a 

In order to identify possible instabilities we select out exponential plane-wave 
solutions of the perturbation equations: SQ = Qae Tt+lkx , where Qo is constant and 
t and x are two orthogonal coordinates in Minkowski spacetime. As our equilibrium 
background state is a fluid at rest we put u a d a = dt- 
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With these assumptions the set of perturbation equations follow as 

= TSn + ikn5u x , 

= TSe + (e + p)ik8u x + ik5q x , 

= T{e + p)5u x + ik(d e pSe + d nP Sn)+TSq x + ikSU xx , 

= T{e + p)5u y + TSq y +ikSU xy , 

= T(e+p)6u z +TSq z +ik5U xz , 



(30) 



Sq x + ik\(d e TSe + d n T5n) + XTTSu x + X-TSq x , 

Sq y + XTTSu y + X-T6q y , 

e 



Sq z + XTTSu z + 

6U XX + ikf]5u x , 
5W V + ikr]5u y , 
SU XZ + ikri8u z , 
5U yy + ikf] v 5u y , 
SIl zz + ikrj v 5u z 7 
5U zy . 



T 

X-T5q z 

e 



Here we have introduced a shortened notation for fj = r] v + |?y. We can put the 
equations above into the following matrix form 

(31) M A B SQ B = 0. 

Here 5Q B represents the list of fields which describe the perturbation of the fluid: 

SQ = (Sn, Se, 6u x , Sq x , SU XX , Su y , Sq y , 5U xy , SU yy 
Su z ,8q z ,SU xz ,8U zz , 8U yz ). 

Then the 14x14 matrix M can be written in the block diagonal form 

/N 0\ 

R 

R 

\0 1/ 

where the submatriccs R and N are defined as follows 



(32) 



M = 



(33) 



(34) 



N 







e+p)T 


r 


ik 












XTT 


i + Arf 












R = 




ikrj 





1 











K 


ikfj v 














/ r 







ikn 









o\ 







r 


ik(e +p) 




ik 







ikd n p 




ikd e p 


T(e+p) 




r 




ik 


ikXd n T 


ikXd e T 


XTT 


1 


f AT? 





\ o 







ikf] 









1/ 
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Exponentially growing plane-wave solutions of (|31|) emerge whenever T and k 
satisfy the dispersion relation 

(35) det M = (det N) (det R) 2 = 

with a positive real T. The roots of this equation are the roots obtained by setting 
the determinants of either N or R to zero. 
The determinant of R gives the condition 

XT-T 2 + ( e + p + k 2 ijX- ) T + fjk 2 = 0. 



The real parts of the roots of this polinomial are negative because the coefficients 
of both the linear and the quadratic term are positive. 

The determinant of N gives the following dispersion relation 

Ap^r 4 +(e+p + k 2 v\^j T 3 + k 2 U + X^(nd nP + P d eP ) - Xnd n T\ T 2 + 

k 2 ((e + p)d e p + nd nP + k 2 fjXd e T) V + k 4 Xn {d e Td n p - d e pd n T) = 0. 

According to the Routh-Hurwitz criteria [3T], the real parts of the roots of a 
fourth order polynomial clqx a + a\x 3 + a^x 2 + a 3 x + a± = are negative whenever 

a > 0, 

ai > 0, 

a\a 2 - a a 3 > 0, 

(36) (0102 — 0003)02 — a±a\ > 0. 



We can see, that the first two conditions of (|36|) are fulfilled according to the 
Second Law, the nonnegativity of the entropy production. 
Let us recall the conditions of thermodynamic stability 

(37) d e T > 0, 

(38) d n ^ > 0, 

(39) d e Td n ^-d n Td e ^>0, 
and the following useful identities 

Td eP = {e+p)d e T + nT 2 d e ^, 

T8 n p = {e+p)d n T + nT 2 d n ^. 

Now, the third condition can be written in a simplified form as 

aia 2 — a a 3 = k 2 X{nT) 2 d n ^- + k^Xif h k 2 f](e + p) + 

1 e 



(40) fc4 ^ 2 § ( d eTp 2 + 2nd n Tp + T 2 n 2 d n ^j > 0. 



T 

The first three terms in the expression are positive. In the parenthesis of the last 
term we can recognize a second order polynomial of p. The discriminant of that 
polynomial is negative 

(41) Di = (2nd n T) 2 - Ad e Tn 2 T 2 d n ^ = -An 2 (d E Td n ^ - d n Td e ^) < 0, 
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because of the last condition of thermodynamic stability. Therefore the expression 
in the parenthesis is positive for all p. 

Hence the fourth condition of (|36p expands to the following form 

(aid2 — 0003)02 — 04a 2 = 

k\e +p)l (d e T(e + p) 2 + 2nd n T(e + p) + T 2 n 2 d n ^) + 

\k i ^(d n T{e+p)+d n ^nT 2 ) 2 + 

2 

Afc 6 (e(e + p)d e T + d e T(e + pf + 2nd n T(e + p) + T 2 n 2 d n ^j + 
X 2 k 6 ^(p(e + p)d e T + n 2 T 2 d n ^) [(p(e + p)d e T + (e + 2p)d n T) 2 + 
n 2 e(e + 2p){d n T) 2 + n 2 T 2 d n ^ (2p 2 8 e T + n 2 T 2 d n ^ + 2n(e + 2p)d n T^j 

A 2 fc 8 77 3 -9 e T + \ 3 k 8 r] 2 dTe 2 {nd n T + pd e T) 2 > 0. 
e 

In the first and third term we recognize the same polynomial expression of (e+p) 
as in (P0]) for p. Therefore all terms are clearly positive, only the term in the 
rectangular parenthesis requires separate investigation. We may recognize that it 
is a second order polynomial of e, with the discriminant 

D = -4n 2 (nd„T + pd e T) 2 (2p 2 (d e Td n ^T 2 - (d n T) 2 ) - (pd n T + nT 2 d n ^) 2 ) < 0. 

The coefficient of the e 2 term is d e Tp 2 + 2nd n Tp + T 2 n 2 d n ^ > is a positive 
quantity according to (|¥T|) . Therefore the term in the rectangular parenthesis is 
positive, too. 

We conclude that the homogeneous equilibrium of the relativistic heat conduct- 
ing viscous relativistic fluids is stable in contrast to the corresponding equations 
of an Eckart fluid. We did not need to exploit any special additional stability 
conditions beyond the well known thermodynamic inequalities and the stability 
conditions of fluids. This is in strong contrast to the Israel-Stewart theory, where 
one should assume additional conditions [15] . 

6. Conclusion 

In this paper we addressed causality and stability, the two most important con- 
ceptual issues in relativistic hydrodynamics. 

We have collected arguments that in dissipativc first order relativistic fluids ap- 
plied to heavy ion collisions acausality related problems may be beyond the range of 
validity of the theory Therefore first order theories, parabolic and mixed parabolic 
hyperbolic equations can be useful physical models in relativistic hydrodynamics. 

Moreover, with a proper distinction of the total and internal energies we sug- 
gested a simple modification of the Eckart theory and we proved that its homoge- 
neous equilibrium is stable by linear perturbations in case of the Eckart choice of 
the velocity field. With the Landau-Lifshitz form of the velocity field our theory 
simplifies to the Eckart theory, therefore the corresponding homogeneous equilib- 
rium is stable, too |22j . However, the Landau-Lifshitz convention has some other 
undesirable properties, that we intend to discuss in a consequent paper. 
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The suggested relativistic form of the internal energy depends on the momentum 
density, therefore the entropy function is the function of the momentum density as 
well. Moreover, in the first approximation we have a regular second order theory 
with only one additional quadratic term in the entropy four vector. However, there 
was no need to introduce additional parameters, the coefficient of the quadratic 
term, and therefore the relaxation time in the generalized Fourier equation, is fixed. 
It was very important, that the general entropy vector is not a simple quadratic 
function in the heat flux. We have got a correction to the viscous bulk pressure, 
too. 
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